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Abstract 

The aim of this paper is to introduce a new technique for calculation of observables, in particular 
multiplicity distributions, in various statistical ensembles at finite volume. The method is based on 
Fourier analysis of the grand canonical partition function. Taylor expansion of the generating func- 
tion is used to separate contributions to the partition function in their power in volume. We employ 
Laplace's asymptotic expansion to show that any equilibrium distribution of multiplicity, charge, en- 
ergy, etc. tends to a multivariate normal distribution in the thermodynamic limit. Gram-Charlier 
expansion allows additionally for calculation of finite volume corrections. Analytical formulas are pre- 
sented for inclusion of resonance decay and finite acceptance effects directly into the system partition 
function. This paper consolidates and extends previously published results of current investigation 
into properties of statistical ensembles. 
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I. INTRODUCTION 

During the last couple of decades the statistical hadronization model, first introduced by 
Fermi [l| and Hagedorn j^, has been surprisingly successful in describing fundamental proper- 
ties of systems created in heavy ion collisions, cosmic rays, and elementary particle reactions. 
In the context of heavy ion collisions this model has been applied to an extensive set of data 



on hadron production, ranging form the center of mass energies o 



the experiments at the SIS, 



AGS, SPS, and most recently, RHIC facilities. For reviews see {3, 4], 0, A systematic evo- 
ution of thermodynamic parameters as collision energy (and size of colliding ions) is changed 
^1 has allowed to establish the 'chemical freeze-out line', which is now a vital part in our 
understanding of the phase diagram of strongly interacting matter. More controversially this 
model has also been applied to a range of elementary collisions Jsj] , where only few particles are 
produced, and the picture of a gas of hadrons can hardly be suitable. The remarkable ability 
of the statistical model to explain these data has lead to the suggestion j^, 9] that thermal 
particle production is a general property of the hadronization process itself, rather than the 
result of a long sequence of microscopic interactions. In this work we will not argue about 
possible physical interpretations [h]] of the partition function of statistical mechanics. However 
we note that in order to apply a semi-classical approximation a volume of O(10fm 3 ) seems to 
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be sufficient 

One of the answers still outstanding in high energy physics is the one of a possible forma- 
tion of a deconfmed state of matter, where degrees of freedom are quarks and gluons, rather 
than hadrons, and the nature of the phase transition between these two phases. The growing 
interest in the study of event-by-event fluctuations in strong interactions is thus motivated by 
expectations of anomalies in the vicinity of the onset of deconfmement 12] and in the case when 
the expanding system goes through the transition line between quark-gluon plasma and hadron 



gas 



13], In particular, a critical point of strongly interacting matter may be accompanied by 
a characteristic power-law pattern in fluctuations [14]]. Multiplicity and charge fluctuations 
have been indeed proposed to be a good discriminating tool between quark-gluon plasma and 
hadron gas 15], provided the signal survives the phase transition and subsequent evolution of 
the system. However, in order to properly assess the discriminating power of such observables, 
one should firstly calculate fluctuations in a hadron gas by including all known physical effects, 
such as conservation laws, quantum statistics, resonance decays, kinematical cuts, finite spatial 
extension, etc. 

Only recently, due to a rapid development of experimental techniques, first measurements of 
fluctuations of particle multiplicities 16] and transverse momenta 17] were performed. And in 
fact one is tempted to interpret recent NA49 data on multiplicity fluctuations in most central 



Pb-Pb collisions [18] as the first observation [19] of the recently predicted canonical suppression 
of fluctuations. The most promising region in the phase diagram for observation of critical 
phenomena seems to be accessible to the SPS accelerator {^(J. A new SPS scan program [21I 
for different ion sizes as well as center-of-mass energies has been proposed to study strongly 
interacting systems at different energy densities and life times. This should be as well our 
main motivation for further investigation of properties of statistical ensembles. The aim is the 
calculation of 'base-line' fluctuations on top of which one hopes to find unambiguous signals of a 
ohase transition [ljj, a critical point [14J], or thermal/chemical (local or global) non-equilibrium 
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The main subject of the past study has been the mean multiplicity of produced hadrons. 
However, there is a qualitative difference in the properties of the mean multiplicity and the 
scaled variance of multiplicity distributions in statistical models. In the case of the mean 
multiplicity results obtained in the grand canonical ensemble (GCE), canonical ensemble (CE), 
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and microcanonical ensemble (MCE) approach each other in the large volume limit. One refers 
lere to as the thermodynamical equivalence of statistical ensembles. It was recently found 
23I . I24I ] that corresponding results for the scaled variance are different in different ensembles, 

and thus the choice of ensemble remains a crucial one, even for large systems. 



In previous publications 
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27 



28 



30t ] the focus was on asymptotic values 



of the scaled variance, or calculations were altogether only performed in the GCE [3JJ. In this 
paper we try to present not only approximations to CE and MCE distributions in the large 
volume limit, but also to find a reasonable approximation scheme for finite system size in an 
analytical, rather than Monte Carlo [32j, approach. The main subject of this work will be 
to identify the GCE partition function with the characteristic function of a statistical system. 
This approach is more effective in terms of mathematical ease as well as in terms of computing 
time than previous methods. Thus we are presenting for the first time analytical formulas for 
the final state of a general multi-specie CE or MCE hadron gas at finite volume. 

The paper is organized as follows. In Section [Til the generalized partition function will 
be introduced and some basic assumptions of this paper will be stated. Section II III gives a 
mathematical approximation to the system partition function in the form of Laplace expansion 
for the CE, while Section [IV] is concerned with the MCE formulation of multiplicity distributions 
in terms of this expansion for a simple quantum gas of massless particles. A general recipe for 
calculation of the scaled variance is presented. Resonance decay is included into the system 
partition function in Section [V] In Section [VI] we show a comprehensive comparison of this 
method to the previously used microscopic correlator approach. A method for finite volume 
corrections is developed in Section IVlII Section IVlIII will present an application of this method 
in great detail and aims to give some physical interpretation. A summary in Section IIXI closes 
the paper. 



II. GENERALIZED GRAND CANONICAL PARTITION FUNCTION 



In textbooks on statistical mechanics (see e.g., Ref. 33J) often first the MCE is introduced, 
where exact conservation laws for energy-momentum and particle number are imposed. Relax- 
ing the constraints for energy and momentum constitutes the CE, while allowing additionally 
particle number to fluctuate about some mean value introduces the GCE. In a relativistic gas 
of hadrons quantum numbers (charges) will be the conserved quantities rather than particle 
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numbers. In this paper it will prove to be of considerable advantage to start off with the GCE 
formulation and imposing exact conservation laws thereafter. The basic idea is to define the 
probability of a given number of particles of some species, iVj, at fixed value of conserved charge, 
Q, i.e. the CE distribution P(N t \Q), in terms of the GCE distributions, P{N h Q) and P{Q). 
Thus, the GCE partition function will be the basis for all calculations in this work. Generally 
the (micro) canonical partition function is obtained from the grand canonical one by multipli- 
cation with Kronecker (or Dirac) delta-functions which pick out all microstates consistent with 
a particular conservation law. It is often more economical to use the Fourier representations 
delta- functions, rather than the delta-function themselves. A short example will motivate the 
following general treatment. The GCE and the CE partition functions Z and Z® are connected 
as, 

Z(V,T,n Q ) = exp(^) Z*(V,T), (1) 

Q=-oo ^ ' 

where V, T, and hq are respectively the system volume, temperature and chemical potential 
associated with conserved charge Q. The probability of finding the GCE system in a particular 
charge state Q equals to the number of all states with net-charge Q divided by all accessible 
states: 

( Q^q \ yQ 

all states with charge Q e ^ \ t ) ^ 
F[Q) = all states = Z ' (2) 

It is important to note that the GCE partition function Z depends on /jLq/T, while the usual CE 

partition function does not. We introduce the generalized grand canonical partition function 

(GGPF) Z{(f>o) by substitution Hq/T — > /iq/T + i(f)Q in Z. One then finds: 



P(Q) = I j ^ e-** Z{<t> Q ) . (3) 

— TT 

The system with one net-charge, Q, could be a quite general one and may include different 
particle species. In the state of chemical equilibrium the chemical potential of any specie / 
equals to /i/ = q\ /j,q, where q\ is the charge of a particle of specie /. Considering only the 
distribution of particle species /, the joint probability P(iVj, Q) to find particle number Ni and 

net-charge Q in the GCE equals to: 

all states with N t particles and charge Q 



P{N h Q) 



all states 

i j d£ j ^ e . iWe e _ imni z(Mni) (4) 
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Similar to Eq. the distribution (@J is presented in terms of the GGPF Z{4>Qi4>Ni) through 
substitutions fii/T — > Hi/T + i4>N t , and Hq/T — > ^q/T + i4>q in the GCE partition function Z. 
Finally the canonical (or conditional) multiplicity distribution P(Ni\Q) is given by: 

all states with Ni particles and charge Q P{Ni, Q) 
^ ^ = all states with charge Q = P{Q) ' ( ' 

Eq. (jSJ) presents the CE particle number distribution P(Ni\Q) in terms of two GCE distribu- 
tions, P(Ni,Q) and P(Q). A detailed account of this is given in terms of an ideal Boltzmann 
pion gas in Section IVIIII 

I I 

It is worth noting that Eqs. (121141151) are as well the basis for any Monte Carlo approach [321 ] . 
A sampling distribution, usually taken from a Boltzmann GCE system, is used to generate 
a {iVj}-tuple of particle multiplicities of all considered species. All 'events' consistent with 
certain constraints, like a set of conserved charges, are accepted, while the rest is rejected. On 
the basis of this set of all accepted 'events' one constructs an ensemble by using a suitable 
re-weighting scheme to account for quantum statistics and proper normalization. One is now 
ready to calculate distributions P(Ni\Q), i.e. the conditional distributions to find particle 
multiplicity Ni in the 'Monte Carlo ensemble', while global charge is fixed to Q, and therefore 
observables like mean multiplicity and multiplicity fluctuations. The advantage of the method 
presented here is certainly that we proceed in a completely analytical fashion and therefore 
much unneeded information, like the exact composition of {Ni}, is 'integrated away'. 

An immediate consequence of Eqs. fl2]H]l5]) is that temperature and chemical potentials appear 
in our formulation of (micro) canonical distributions (as well as in the Monte Carlo approach 



321] ). At first sight this seems to be a serious problem and an unnecessary complication of 
our initial task to find a reasonable approximation to the CE and MCE partition functions. 
However, the main technical challenge when numerically integrating the original version of 
the microcanonical partition function arises from a heavily oscillating integrant. Auxiliary 
parameters T and /i will produce a very smooth function, for which approximation schemes 
can be used (see Section [Villi) . In fact chemical potentials and temperature can be factored 
out, and thus our partition function is just the original partition function times some factor 
(Appendices [B] and [C]). Taking the ratio in Eq.flHD artificially introduced temperature and 
chemical potential drop out. The quality of our approximation on the other hand will crucially 
depend on their choice. We will show in Section IVIII that the requirement of maximizing the 
GGPF at some given equilibrium point leads to a unique determination of thermal parameters 



and moreover constitutes the optimal choice for our approximation scheme. This prescription 
is self-consistent and can be shown to be not in contradiction to basic thermodynamic relations 
known from textbooks |33|. 

Taking the GCE partition function as a basis, we will throughout this paper use the lan- 
guage of statistical mathematics for calculation of conditional multiplicity distributions in any 
ensemble. So it is worthwhile to give a general outline. For any probability distribution function 
(PDF) one can define an associated characteristic function (CF) by Fourier back-transformation 
(Appendix In this work we will identify the GGPF as the CF of any GCE distribution, and 
find joint GCE distributions, e.g. P(Ni, Q), by Fourier analysis of the GGPF. In the example 
above the CF associated with the PDF P(Ni, Q) in Eq.fll]) is Z(4>q, tfix^/Z, while the CF asso- 
ciated with the PDF P(Q) in Eq.(j2]) is Z(<J)q)/Z. However, only in the simplest cases one can 
find an analytical solution for these Fourier integrals. It seems as well to be rather difficult to 
obtain the CF directly for the PDF P(N t \Q), and hence we define P(N t \Q) = P(N h Q)/P(Q). 
Provided one has knowledge of all (sufficiently many) moments of an PDF one can construct 
the CF and find the exact (approximate) PDF. For practical purposes, however, cumulants, 
which are obtained by Taylor expansion of the logarithm of the CF, are far better suited. 

III. CANONICAL ENSEMBLE 

The method presented here will provide the basis for the following discussion. We present 
in detail our calculation of the CE partition function. In this section we will employ saddle 
point expansion to approximate the canonical PDF for an ideal hadron resonance gas and focus 
on the asymptotic solution only. We first calculate the 3-dim distribution of charges in GCE, 
P(Q, B, S), namely baryon number, strangeness and electric charge. In fact our calculation 
is similar to that of Becattini et.al. [29], however we start explicitly from a GCE partition 
function which includes chemical potentials. We will identify our version of the CE partition 
function as the un-normalized probability distribution function of conserved charges in GCE, 
and show that this type of integral over the GCE partition function generally leads to a Gaussian 
(multivariate normal) distribution in the large volume limit. In a next step we use this result to 
obtain the 4-dim PDF for GCE, P(N, Q, B, S), hence the probability to find our GCE hadron- 
resonance gas in a state (N,Q, B, S). The ratio is the canonical PDF P(N\Q, B, S), see Eq. 
(jSJ). In particular we will present a very simple formula for the asymptotic scaled variance of 
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multiplicity fluctuations. 



A. CE Partition Function 



The canonical system partition function of a hadron resonance gas with three conserved 
charges, i.e. electric charge Q, baryon number B, and strangeness S is given by the triple 
Fourier integral over GCE partition function {34} ]: 



Z Q. 



B,S 



2tt 



>B 



2tt 



2tt 



e -iB4 Be -iSto x exp 



E 



2/ 



>B, 9S) 



I 



(6) 



Adopting vector notation for the set of conserved charges (Q, B, S) = Q = Q 3 ' and angles 
<ps) = 4> = <Pj we can write Eq.Q in compact form: 

3 



z Q3 



n 



W) 



-iQ 3 



exp 



E 



2/ 



(7) 



where repeated upper and lower indexes j imply summation over j. The single particle partition 
function of particle specie I is given by: 

9 iV r _ / n ±1 



2; 



(27T)" 



y d 3 p In ± e - (£! ~ w)/T e^)" 



V^, (0. 



where we have introduced particle Vs quantum number configuration qj 



3) 1 



(li 



(8) 



(qi,k,si) 



it's degeneracy factor g\ = (2Ji + 1), internal angular momentum Ji, mass mi, and energy 



si = yjp 2 + mf, the chemical potential vector fi 3 = (//q, /j,b, fJ-s), and particle Vs chemical 
potential fii = qj\iy V is the system volume, and T it's temperature. The summation 
includes also anti-particles, for which qj — > —qj . The (Wick rotated) fugacity of particle specie 
I is thus given by the substitution \ t = exp [qjfij/T] — > \i = exp \qj (f-ij/T + i</>j)]. The upper 
sign in the Eq. (jSJ) denotes Fermi-Dirac statistics (FD), while the lower is used for Bose-Einstein 
statistics (BE). The case of Maxwell-Boltzmann statistics (MB) is analogous and will not be 



discussed separately. Together with Y^i^i 
written as: 



the system partition function can be 



Z Q1 



n 

3=\ 



LUfJj 

W) 



exp 



V * 



(9) 



9 



The following calculation is based only on the general form of Eq.flHJ), and can easily be gen- 
eralized to a larger set of conserved quantities. In the language of statistical mathematics the 
function ^ (<f) j) would be called cumulant generating function (CGF) of the partition function 
Z® 3 , while exp [V ^ ((f) j)] is called CF. For large volume, V — > oo, the main contribution to the 



integral Eq.Q comes from a small region around the origin 



291 ] . Thus we proceed by Taylor 



expansion of the CGF \& ((f) j) around <j>j = and introduce cumulant tensors k: 



31 ,32, — ,]n 



-I 



(10) 



The CGF therefore can be expressed in terms of a Taylor series: 



n=0 



J2 • • ■ rjn J 



where summation over repeated indices is implied. The cumulant of th order is just the 
logarithm of the GCE partition function Z divided by the volume, Z = exp(V/t ). Hence, after 
extending the limits of integration to ±oo, which will introduce a negligible error, we find: 

3 



z QJ 



n 

3=1. 



W) 



exp 



— i 



Q j 4>j + vY,-<' j2 '-' jn ^h 



71=1 



y J2 • • • VJM 



(12) 



It is worth noting that the integrant of Eq. (|T2|) is not 27r-periodic anymore, while the one of 
Eq.Q is. Spelling out the first two terms of the summation yields: 



z Q] 



n 

3=1. 



(2 



7T 



exp 



'3\vn 



32, — ,]n 



n=3 



(13) 



Performing now a change of variables will simplify this triple integral. 



V a/ (f) k , 



where <7 fc is the square root of the second rank tensor k,2'- 



(7 \ 



The element ddj equals to: 



detlv^a] d(f)j = V 3/2 det\a\ 



(14) 



(15) 



(16) 
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Lastly in terms of this transformation normalized cumulant tensors A are introduced: 



pi, j2 ,..,^ = K 



kl,k2,—,k, 



(17) 



The new variable £ 3 ' will be a measure for the distance of the actual charge vector Q J to the 
peak of the distribution of the PDF: 



Including all these steps at once yields: 
Z 



Z Q3 



V 3 ' 2 det I a I 



3 

n 

3=\ 



dB 3 
27 



exp 



WO, 



2! 

\ixj2,."jV, 
.n +1 A n 



n=3 



r/: 



^il %2 • • • @jn 



(19) 



Eq. (TT9T) is the starting point for obtaining an asymptotic solution in this section as well as for 
finite volume corrections in Section IVHI Through coordinate transformation Eq. (1141) we have 
explicitly separated terms in their power in volume. Thus as system size is increased influence 
of higher order normalized cumulants A ra decreases, allowing for truncation of the summation 
for sufficiently large volume. A few words on physical units are in order. The single particle 
partition function Eq.® tpi[fm~ 3 ], and therefore all cumulant elements Eq. ffTO"]) «;[/m~ 3 ] in 
CE, consequently entries in Eq. (ll5l) a[fm~ 3 / 2 ]. The normalization in Eq. (fl9l) V 3 ^ 2 det |<r| for 
3-dim a, as well as the new variable of integration Eq. (1141) J , are hence dimensionless. The 
inverse sigma tensor elements are a~ l [fm 3 / 2 ] and thus the charge vector Eq. (jl8p £ J will be 
dimensionless. Finally normalized cumulants Eq. (1171) are X n [fm~ 3+3n ^ 2 ], which is canceled 
by the factor V~ n ^ 2+1 in the summation in Eq. ffT9l) . Thus all terms involved in Eq. ffT9~]) are 
dimensionless. 

For V — ► oo one can discard terms of V~ x l 2 and higher: 



Z Q3 



Z 



V 3 ' 2 det Id 



3=1 



exp 



6% 
2! 



(20) 



Completing the square, the integral (1201) can be solved: 

5 J tj 



exp 



{2nV) 3/2 det \a\ 



(21) 
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From Section [TT1 it follows that in GCE one would find for the 3-dim charge PDF: 



P(Q ] ) 

Please note that Eq. (l2"2"j) (a 
scopic correlator approach 



e t 



cxp 



{2nVf /2 det\a\ 



(22) 



beit in different notation) was used as an assumption in the micro- 



24 



25 



2a,|23,|28 1 . 



B. Particle Number Distribution 

Similarly to Eq.Q the canonical state involving a charge vector Q- 7 and a particular number 
of particles N of some species can be described by the following 4-dim partition function: 



Z n,Q,b,s = / am / u*PQ / u<PB / ays iN4>N iQ(t>Q m <t, B -iS<t> s 
2tt / 2vr / 2vr / 2vr 



x exp 



(<f) N , (f)Q, (j) B , 0s) 

.1=1 



(23) 



where (N,Q, B, S) = (N,Qi) = Qi . Each particle specie which is selected receives a Wick 
rotated fugacity exp («0at) in complete analogy to those associated with charge conservation. 
Similar to Eq. fl21l) the canonical partition function Z® 3 can be approximated by: 



exp ' 



Z Q3 ~ Z \r n '— , (24) 

(2vr\/) 4/2 det|a| 1 ' 



Therefore in GCE we find the 4-dim PDF of finding our system in a state Qi\ 

ii c. 

2?Q> exp ' - " 



However we are only interested in a 1-dim slice with constant values of conserved charges Q J = 
(Q,B,S), e.g. the conditional PDF P(N\Q j ) = Z& /Z Qj , where the normalization is Eq.(l2H 
from the previous section. The expansion works best around the peak of the distribution. If 
analytical solutions were available this would not matter much, but here we need to find the peak 
of the PDF. In the thermodynamic limit peak and mean of the distribution coincide (see Section 
IVTT]1 . hence it peaks at Q k eq = Vk\ = ((Q), (B), (S)). For the 3-dim distribution Eq.(l2Tj) from 



the previous section it follows that Eq.(Tl8]) £ j = (Q k eq - Vk\) (a 1 ) k j V 1/2 = 0. For the 4-dim 
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distribution of this section we find at Q k eq = (N, Q k eq ) for = (N — (N)) 2 V~ l (a" 1 ) / (o^ 1 )/, 
and (N) = Vnf. The canonical PDF can thus be written as: 

Hence one can find the width of the distribution from its normalization, as well as from its 



exponential det |a|/det \a\ = \/ {a' 1 )^ (a- 1 )/ = ((N 2 ) - (iV) 2 ) -1 / 2 V 1 ' 2 . The identity 

(^/r 1 ); ■ (27) 

is proven in Appendix ([D]). One gets for the scaled variance of the particle number distribution: 

(N 2 )-(N) 2 V det\a\ 2 det\k 2 \ 



det 




det 





In words, this is the ratio of the product of the eigenvalues of the 3 + 1 dimensional matrix Ri and 
the 3 dimensional matrix k 2 divided by the particle density of the particle under investigation. 
The equivalence of Eq. (l28l) to results of the micro-correlator approach, see e.g. 24, 25], 



26 



27 



28] , will be shown in Section I VI I 



\2 



wSjp-'i-W 1 • (29) 

Thus we obtained the rigorous mathematical proof of the hypothesis that the multiplicity dis- 
tributions in a hadron-resonance gas in thermodynamic limit have Gaussian shape. Diagonal- 
ization of k 2 was not at all necessary for arriving at the asymptotic solution Eq. (f28|) . However 
for finite volume corrections one will have to go through the process of finding eigenvectors and 
eigenvalues of k 2 (see Section |VTT]). 



C. The Cumulant Tensor 

In this section we will calculate the K-tensor Eq. fflO|) only for primordial single particle specie 
fluctuations. In Section |V] this will be extended to fluctuations of a selection of particles, like 
'positively charged' and resonance decay will be included. 



= f -i J- \ * 



(30) 

^=0 



13 



and 



K 



31,32 



d 



y 3l 



d 



—I 



32 



4>i=0 



(31) 



Cumulants of order 1 give GCE expectation values, hence average baryon, strangeness, electric 
charge, and particle density. The second cumulant contains information about GCE fluctuations 
of some quantity (diagonal elements), as well as correlations between different quantities (off- 
diagonal elements). For the distributions P(Q j ) and P(Q j ) the first and second cumulants 
are: 





f K Q,Q 


Q,B 


K Q,s\ 


K 2 = 


B,Q 
K 2 


B,B 
K 2 ' 


B,S 
K 2 > 




S,Q 
\K2 


S,B 
K 2 


S,S 
«2 / 



«2 



= K 


Q 
K 1 , 








iV,Q 
K 2 


N,B 
K 2 


N,S\ 


hu 2 


hu 2 


K Q,B 

r"2 


hu 2 


B,N 
K 2 


B,Q 


B,B 
K 2 


B,S 

hu 2 


S,N 
\«2 


hu 2 


S,B 
K 2 


s,s 

K 2 1 



(32) 



(33) 



In case two quantities are un-correlated, as for example primordial 7r + multiplicity and globally 
conserved strangeness (n + does not carry strangeness), then the corresponding elements k 2 ,s = 



S N 

k 2 =0. For clarity some elements are explicitly given. The primordial mean value of particle 



number density of particle species / and the mean charge density are: 

d 



N 



—I 



—I 



d<j>N 
d 



3 =O 



3 =O 



(34) 
(35) 



Fluctuation of particle density k^'^, correlation between particle number and baryonic charge 



and correlation between strangeness and baryonic charge k 2 ' B are given by: 



N,N 

t\,r) 



N,B 

tvn — 



—i 



9 Y* 








d(fi N J 


3 =O 


9 W 







—i 



_d_ 



d 



—i 



- = 6i < , 
=o 



where the first and second derivative of ip are: 

91 i 



< = 



(2tt)< 
(2tt) j 



d 3 p 



e -(e,- w )/T 
(1 ± e -(e«-w)/T) ' 
e -(ei-m)/T 

(1 ± e-fe-«)/T) 2 • 



(36) 
(37) 
(38) 

(39) 
(40) 
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Please note that the analogs of the matrix K2 and the vector k%, however in different notation, 



were used in both previously published methods for calculation o 



24 



scaled variance under the 



25 



26 



21, 



281 ] . and saddle 



thermodynamic limit, the micro-correlator approach, see e.g. 
point expansion method . The advantage of this paper is certainly the very simple formula 
(128)) . and the possibility for calculation of finite volume corrections (see Section l"VHI) . giving a 
lower bound for the validity of the micro-correlator approach and saddle point expansion. Exact 



agreement was 



limit 



19 



24 



25 



bund with the analytical results for the scaled variance in the thermodynamic 



2fl 



27 



281 ] . Numerical calculations give very good agreement to [29|, despite 



the fact that different particle tables were used in the latter. 



IV. MICROCANONICAL ENSEMBLE 



In the MCE we additionally enforce kinematic conservation laws. The simplest example is 
an ultra-relativistic gas, made up of one kind of neutral massless particles. The reason for this 
choice is that an analytical solution exists, at least for a microcanonical ensemble, where energy 
but not momentum is conserved [l[ [27], in MB approximation, allowing for a comparison of 
analytic results to asymptotic solutions (see Sections IVH CI and IVHI Cl) . Further we can carry 
out the integration of the GGPF in full detail for this example (Appendix [C]), and highlight 
some features of our method, such as the role of temperature in the MCE. 

In this section, however, we explicitly include momentum conservation as well as FD and 
BE statistics, and restrict ourself to the asymptotic solution. The basic ideas are essentially the 
same as before. The only conserved quantities here then are the total energy E and the three 
momentum P. In principle we would have to treat energy and momentum conservation on equal 
footing and introduce Lagrange multipliers associated with conserved momenta. However, in 
the rest frame of a static thermal source we always find (P) = 0. Thus we only need to consider 
the th component, 1/T, of the 'four-temperature' [32]. The probability to find a microcanonical 
system in a state with exactly iV particles is defined by: 



P{N\E,P = 0) 



number of all states with N particles, E, and P = 



number of all states with E, and P = 
Our starting point is again the GCE partition function, which can be written as: 

9 



(41) 



exp 



V- 



(2tt)< 



d A P f (p) 



(42) 
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where f(p) is the probability of having a particular momentum state occupied. For massless 

p |= a/Px + Py + Pzi thus f(p) = exp (— \p\/T) in MB approximation, 



particles we find e 

or in the quantum statistical treatment f(p) = m(l ± exp (— \p\/T)) 1 , where the upper sign 
denotes FD, while the lower sign stands for BE statistics. The numerator in Eq. (1411) is given 
by the 5-dim Fourier integral over the GGPF Z i y> N , ip E , y> p 



PIX.E.P) =i 



oo oo 

ijv f d<t>E 



2tt 



2tt J (2 



^ e -W^ e -^fe e -iP4> P exp 



^7T 



'AT, <PE, <Pp 



(43) 

Energy and momentum conservation (continuous quantities) require additionally the use of 
Dirac 5 functions, rather than the Kronecker 5 for discrete quantities. We need dimensionless 
quantities for the delta-functions, thus 4> E and </> p need to have dimension [Ge^ _1 ] (see as well 
Section IVlII C\ Appendix (Cj and ref. |32j). The CGF can be expressed with a suitable choice 
for cf>E and <j) p in MB statistics as: 



>JV, <PE, 9p 



(2*)' 



d 3 p e~^l T e i( ^ N e i '^^- B e^ p 



(44) 



or for FD (upper sign) and BE (lower sign) statistics as: 



>N, (PE, 9p 



(2^ 



d 3 p In 



±i 



(45) 



Using shorthand Qi = (N, E, P x , P y , P z ) and <j>j = [<f)jsr, 4>Ei <fip x , 4>p y i 0p z ) one can simplify the 
notation, and only use the asymptotic solution Eq.(l22l) for large volumes, derived in the previous 
section. 



P Q 



Hj (27T) 

J — CY1 



1 exp 



e -iQ 3 4>j 



ex p \ V J2 % -\ <"' ]n ^ 



n=l 



(46) 



(2nVf /2 det\a\ 

The cumulants needed are now given by the respective derivatives at the origin. In case FD or 
BE statistics are used one will employ the same set of derivative operators, which result in the 
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usual, yet slightly more difficult, integrals due to the logarithm. The expectation values are: 



.N 







—l- 







—l- 



d 



J Px 



<f>i=0 



4> j= o 



(27T)' 

9 
9 

(27T) 



<i 3 p /' (p) bl , 



d 3 p f (p) = 



(47) 
(48) 
(49) 



while some selected elements from the second rank tensor K2 are: 

2 

'<Ppf"{p) 



N,N 



d 



5 



K 



„Px,Px 



E,N 
K 2 







a 



—i- 



d 



4>j=o (2-kY 
9 

4>i=S ' (2tt) 3 
9 

* s =8 (2tt) 3 
d 



d 3 pf(p) |pf 



d 3 p f" (p) p£ 



where similar to Eqs. (139]) and fj40|) we define: 



(2tt) 



d 3 P r (p) \p\ 



(50) 
(51) 
(52) 
(53) 



f(p) 



-\P\/T 



and 



/" (p) 



-Ip1/t 



\2 ' 



(54) 



1 ± e-|pl/ T ' ' vw (l ± e-|pl/ T )' 

or simply f {p) = f (p) = f" (p) in MB approximation. Due to spherical symmetry of the 



momentum distribution, we find K2 X,Px 



J>v>Py 



K 



PZ;PZ 



. Further correlation terms are identical 



to zero, k 2 



N,p x 



E,Px 



JPx,Py 



0. In our previous notation we find the first and second 



cumulants, k\ and k 2 , in Boltzmann approximation for the distribution P(N, E, P): 

(snT 3 , 24vrT 4 , 0, 0, o) , and (55) 



«2 



(2^ 



( 8ttT 3 24ttT 4 
24ttT 4 96ttT 5 
32nT 5 









\ 















32?rT 5 

Q n QO^T5 



(56) 



327rT 5 / 



For the normalization we need additionally the second cumulant k 2 of the PDF P(E, P), which 
is obtained from R2 by crossing out the 1 st column and I s * row. The remaining relevant integrals 
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MB 


FD 


BE 


< 

N,N 
K 2 

E,N 
K 2 

E,E 
K 2 


8vrT 3 
24vrT 4 

8vrT 3 
24vrT 4 
96vrT 5 


6ttC(3)T 3 
18vrC(3)T 4 

14 5 T 5 
15" 1 


8ttC(3)T 3 

±^ 

§vr 3 T 3 
24vrC(3)T 4 

16 5 T 5 
15 /l ^ 



TABLE I: Selected elements of the first two cumulant tensors for a massless gas in Fermi-Dirac, Bose- 
Einstein statistics and Boltzmann approximation. All entries have to be multiplied by gj (27r) 3 . The 
Riemann Zeta function is C (3) — 1.202. 

for FD and BE statistics are summarized in table [JJ Following the recipe from above, see 
Eq. fT28l . the scaled variance can be expressed as follows: 

2 

(57) 



E,E N,N 



UJ 



det 


«2 




< det | 


K2| 



I Kj 



N E,E 

b 1 hu 2 



Independent of energy density one finds for the large volume limit the following asymptotic 
values for the scaled variance u, in exact agreement with 27|, in table [Til Momentum conser- 





MB 


FD 


BE 


^mce 


0.25 


0.198314 


0.535463 



TABLE II: Asymptotic of the scaled variance in a MCE for a neutral massless gas 



vation explicitly drops out in the thermodynamic limit in the calculation of the ratio of the 
relevant determinants. This will not hold true for any finite system size (see Section IVIIj) . 



due to the appearance of off-diagonal terms symmetric in momentum like k^' Px,Px . However 
for large enough volumes, expected to be created in heavy ion collision experiments, one can 
probably safely disregard exact momentum conservation, when considering multiplicity fluctu- 
ations in the full momentum space, i.e. 47r yield and fluctuations. In previous publications on 
asymptotic multiplicity fluctuations in the MCE it was only always assumed (without proof) 
that exact momentum conservation would not affect the result, and thus for technical reasons 
not taken into account. 
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One last issue should be mentioned. For partition functions for very small volume one 
should return to summation over (quantized) momentum states, rather than integration over 
(continuous) momentum space, or ideally turn to a quantum field theoretical frame work of 
the MCE |ll|]. Usually these distinctions do not affect the result much, simplify however 
calculations considerably For a massless gas however we encounter the problem of a divergent 
k,^ ,N,N for BE statistics. This problem is well known from textbooks and usually overcome by 



applying a low momentum cut-off 33]. 



V. RESONANCE DECAY 

In this section we want to show how to include resonance decay analytically directly into 
the system partition function. This will prove far more efficient than the definition of a gen- 



erating function 26] which requires a rather cumbersome calculation of all possible primordial 
correlators, see Section IVI Bl Particle decay is itself a random process. Nevertheless one can 
assign a particular volume in phase space, given by the value of its single particle partition 
function zi, Eq.([H]), to one type of resonance I. Resonance decay will now populate this volume 
in phase space according to some weight factor, the branching ratio, for each of the possible de- 
cay modes. This weight can be assigned to the particle type(s) we are set to investigate. Based 
on the assumption that detected particles are drawn in the form of a random sample from all 



final state particles, e.g. disregardin g co rre 



acceptance scaling employed in ref. 19 



23 



ation in momentum space, this procedure leads to 



26J (see Section IVl CI) . Conservation laws can be 



imposed on the primordial state (rather than the final state), since decay channels, which are 
experimentally measured, do not only obey charge conservation, but all relevant conservation 
laws (omitting weak decays). 

A. Final State Partition Function 

For calculation of final state distributions, we have to determine final branching ratios of 
a resonance type into only stable particles. As an example we consider the decay channel 
A —>■ B + X with branching ratio r ' a-^b+x — a - Resonance B could itself be unstable and 
subsequently decay via the channel B —>■ Y + Z with branching ratio Tb^y+z = b. So we 



define the final branching ratio Ta-*x+y+z = ^a-^b+x ^b^y+z = a>b. Decay tables in [19L 126] 
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have also been generated according to this prescription. 

For resonances it seems economical to define further absolute branching ratios 1^ as the sum 
over all final decay channels of resonance I with a given number c of selected daughters. Hence 
17 is the sum over all final decay channels with two daughter particles which are of interest, 
i.e. two positively charged particles in case one wants to calculate uj + . As a consequence of 
this definition, branching ratios 1^ will depend on which to one is set to calculate. 

For the final state one has to take all C\ absolute decay channels of resonance type I into a 
number c of selected stable particles into account. All parent resonances carry fictitious particle 
specific fugacity \^ taken to the power of the number c of selected daughters in a particular 
channel, (Xn)°- For the sake of a common treatment for all particles and resonances we 
assign a 'decay' channel to stable particles as well, either = 1 (while T^ 1 = 0) if selected, 
or r° = 1 (while = 0) if not selected. The single particle partition function reads after 
substitution (Xn) c — ► « 



9i 



(2tt)< 



d 3 p\n [ i±e~^-^/ T e iq '^ 



^ ^ pc e ic(f> N 
c=0 



±1 



(58) 



where the form of the vector q\ depends on our choice of ensemble. For instance, we could 
have for a hadron resonance gas with three conserved charges q\ = (qi,bi,Si) in the CE, or 
qj = (qi, bi, si,Ei,pi) in the MCE. The sum over all decay channels which produce from zero up 
to a number Ci particles of the selected types needs to be one: 



c, 

En 

c=0 



(59) 



This is somewhat of a practical challenge, since decay chains of heavier resonances are not 



always well established 



35| and respective thermal models codes |3^, |37|] struggle to implement 



this. There are several ways to deal with this, the two extreme ones are 1) rescale all known 
channels according to Eq. fl59|) . to unity, or 2) assign the missing fraction to the 'channel' T^, 



e.g. to the channel without stable particles of interest. In case one sets all angles 
i.e. one returns to e lc ^ N —* {Xn) c , one obtains the GCE partition function: 

±1 



exp 



9i 



d 3 p In lie' 



-( E! - W )/T 



' Ci 
.c=0 



N, 



(60) 



from which GCE expectation values can be calculated as (N) = A at In Z 



=i- 
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B. Limited Acceptance 



An imperfect detector will modify our definition of branching ratio Tf. Based on the as- 
sumption of un-correlated particle detection, all particles have equal probability q of being 
observed. The corresponding acceptance distribution of observing n particles when c are pro- 
duced is given by a binomial distribution P aC c(n, c) = q n (1 — q) c ~ n r). We define the effective 
branching ratios r™ , which already include the effect of finite acceptance, as the product of 
absolute branching ratio 1^ and acceptance distribution: 



1 1,9 



£ <f (i-g) c 

c=n 



(61) 



For example the effective branching ratio Tf is the sum over all absolute branching ratios T 
which produce at least 2 stable selected particle, Tf q = Tfq 2 + 3Tf q 2 {l — q) + 6rfg 2 (l — q) 2 + . . . 



The binomial coefficient (^) 



takes care of the fact that particles are indistinguishable. 



n!(c— n)\ 

In our approximation resonance decay and particle detection are thus two independent random 
processes. When both, particle decay and un-correlated detection, are included, the single 
particle partition functions reads: 



ipi, g (0j,07v,r^,g) 



9i 



d 3 p In 1 ± e 



(27T) 3 

The sum over all effective channels Tf q is equal to one 

Ci 



E ix. 

n=0 



±1 



(62) 



Ci 

En 



E r ?E' n (i- 9 r"Q = 1 



(63) 



n=0 c=0 n=0 

since the second summation (52n=o) * s ec L ua l to unity, while according to Eq. (l59l) the branching 
ratios 17 are also normalized. The final state CGF thus is the sum over all such single particle 
partition functions: 

* = E M^>^,ry . (64) 



The Cumulant Tensor 



The procedure is essentially no different from Sections II I II and IIV} yet slightly more com- 
plicated, due to the various multiplicities in the decay modes. For clarity some elements are 
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explicitly given. The choice of ensemble then naturally defines the cumulants needed for cal- 
culations. Both, CE and MCE, are considered. The final mean values of particle and energy 
density are: 

c, 



d 



—i- 



<t>j=0 



E 

i 

5>. 



.71=0 



1,9 



77 



(i,i;0) 



(65) 
(66) 



Fluctuation of particle and energy density, and correlations between particle and baryon num- 
ber, particle number and energy, and between baryon number and energy are given by: 



NM 



K ■ n 



E,E 
K 2 



N,B 



d 



>N 



d 



—i- 



—i- 







<t>j=0 

d 



E 

i 

E* 



o. 



Er 

n=0 



n a 
l,q n 



(i,i;0) 



c, 



E r ". 



n 



.71=0 



(2,2;0) 



(1,2;2) 



— Z- 



N,E 



d 



d 



d 



d 



J= O 



fc=0 



' Ci 

E r ;:, 

n=0 



n 



E 



Ci 



E r ". 



n=0 
(1,2;1) 



(1,2;1) 



(1,2;0) 



(67) 

(68) 

(69) 

(70) 
(71) 



where we use a more general shorthand notation for the derivatives similar to Eqs. (l39p . and 
TDD: 

n, f ( e -(ei-Vi)/T) a 



(a,b;c) 



(±1) 



(72) 



(2tt) 3 7 (1 ± e-( e '-w)/r) 6 ' 

It is quite easy to see that higher orders of these elements become quickly very complicated. 
For calculation of asymptotic multiplicity fluctuations, integrated over full momentum space 
(Air fluctuations), it is not necessary to take momentum conservation into account (see Section 
rVT) . The implementation of proper kinematical cuts will be the subject of a forthcoming paper. 



VI. CORRESPONDENCE TO MICROSCOPIC CORRELATOR APPROACH 



One can make a correspondence between different notations in this paper and in the micro- 
correlator approach published previously 
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25 



26 



27. 



28| for asymptotic fluctuations. 
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The aim is to show that both methods lead to the same result, as well as handle resonance 
decay exactly the same way. This fact should not be taken as trivial, as both methods appear 
very differently (not only in terms of notation). In this work all quantities are derived from 
a macroscopic partition function, while in the microscopic correlator approach one considers 
average occupation numbers, and their fluctuations, of individual momentum levels and cor- 
relations between different levels. We restrict the consideration to CE with three conserved 
charges. Additional charges or energy (momentum) conservation in the MCE would just lead 
to enlargement of the matrixes k 2 and k~ 2 . We start off from Eq. fl28|) : 



det 


k 2 




< det 


K2I 



(73) 



where det \k 2 \ and det \k 2 \ are the determinants of the matrices Eq. (l3"3"|) . Expanding the deter- 
minants in terms of their complementary minors, and M^k respectively, yields 38] : 

4 

det |5| = Yi-iy +k Mj,k, (74) 



where minors of k 2 can be expressed in terms of the minors of k 2 : 



M N:N = det |« 2 | , 
After straig! 



M 



N,Q 



Kf N M .o - i£'"M B n + k^M : 



B.N 



S,N 



B,Q 



l S,Q , 



etc. 



(75) 



W c .e. — 



;htforward calculation Eq. (|73l) can be rewritten as follows: 

lr lr ~ at^~ ^rD~ N,S77 

'9 M N)S 



N,NTf N.QTf N,B7~f A 

K 2 Mn^n — «2 -MjV.Q + ^2 Mn,B + K 2 



N,N 
K 2 

~ ~N~ 



det I k 2 I 



o N,Q N,S A/T 

+ 2k 2 '^k 2 ' M QjS 



Ki det I k 2 

K 2 

o N,Q N,B I* 

- 2k 2 * k 2 M QiB 



v m Q:Q + (4> b Ym b , b + 



N,B N,S xyr 

- 2k 2 k 2 M b ,s 



(76) 



This is the most general case. For a specific calculation we need to specify the matrix elements. 



A. Primordial 



In order to change our notation here to that of the microscopic correlator approach JJJ, [26] 
we need the derivatives Eq.( 172l) of the single particle partition function ipf. 

T//1» 1 ;0_V~ S / \ ta/1,2;0_V^ 2 j ;2,2;0 _ / 1.1:0 / 1,2:0 



V 
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where the subscripts i and p denotes particle species and momentum level respectively. The 
matrix entries from Section IHI CI are: 



p. I 



p.'l 



p,l 



Vk$> Q = A(g 2 )^ J2^ v h Vk i' B = E*K* etc 



p. i 



p. i 



where 



33j ] = (n Pi i)(l ± (upper sign for FD, lower sign for BE). It is easy to see 



that our matrix K2 corresponds to the correlation matrix A = VK2, see e.g. Eq.(12) in [26] 
Substituting into Eq.( 176|) one finds Eqs. (50,51) from reference [261 ]: 



UJ, 



with the correlator (A7VjAiVj) c . e , being Eqs. (11, 15) from 261 ]: 



(AJViAAf,), 



71^ ^.2 

v p ,i Oij $pk r^r- 



qtqjMQ^Q + bibjMB,B + SiSjM, 



s,s 



+ (QiSj + QjSi) M Q;S - (qibj + qjbi) Mq B - fas, + bjSi) M B „ 



Eq. f!77|) can be simplified for the case of single-specie fluctuations: 



-r, . = ^g.c.e. 



+ 2g jSj MQ i5 - 2qjbjM Q , B - 2b jSj M BtS 



which coincides with Eq.(16) from 26], where uj 3 g ce = v p J J2 p ( n p,j) 



(77) 



(7S 



(79) 



B. Resonance Decay 



In order to account for resonance decay we need the matrix elements stated in Section IV CI 
together with: 

Ei> 2 = EE< n ^' ^d (Er>) =EE< n *>^>«' ( 8 °) 



n=Q 



* 3 



x n=0 



* 
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where the summations are to be taken over all stable hadrons i and j which one is set to 
consider, and n = rii is the total number of selected particles in this channel. The matrix 
element Vk£' = A (n 2 ) then splits into summations over stable hadron i, j and resonances R , 

A(n 2 ) = J2 v tv + X>L>X>i>*K>i* + J2 nR >pJ2( AniAn ^ R ' ( 81 ) 

i,p R,p i,j R,p i,j 

while the correlation terms Vk 2 = A (nq) are of the form 

A(ng) = 5^g<«Jp + J2qRV RiP ^2(ni) R . (82) 

i,p R,p | 1 i 

Substituting into Eq.(J76l) gives Eq.(47) in reference [26(: 

(ANiANj)^ = (AN*AN*) c . e . + ^(N R ) (An, A nj ) R + ^(AiV* AN R ) c . e . (n,) R 

R R 

+ ^(AiV* AN R ) c . e . {m) R + ^(AiV^ AN R ,) c . e . (m) R ( nj ) R > . (83) 

R R,R> 

Correlation terms (AN*AN R ) c . e ., and (AN R AN R >) c . e , appear in the CE (and the MCE) due to 
products of A (nq) A (ns), etc., in Eq.( |76l) . and are absent in the GCE. 



C. Acceptance Scaling 



Starting off again from our approximation of un-correlated particle detection, from Section 
IV B[ we find for the first two moments 38) of the binomial distribution of detected particles 
produced by decay of resonance type I: 



c, 



)i = J2 T i E nq{l -tf~ n { C n) = E 17 qC = q{c)l 

Ci 

2\ ... -pc „2 „ /-i „\c-n / ° \ _ -pc „ /-i \ i 2 2 



in 



c=0 n=0 



2 )' = E r * E n2 ^ 1 -^" n ( J = E r * [9(l-<7)c + g 2 c : 

c=0 n=0 ^ ' c=0 

= g(l-g)(c), + g 2 (c 2 )/ . 

Comparing with Eqs. (l65|) . (I69I) . and f|67|) we find for the cumulant tensor elements 



K f = E q ( Cj ) ^ 



(i,i;0) 



Q («f)4jr 



and 



B,N 
K 2 



N.N 



E 9 &/ ( q ) ^ 



(1,2;0) 



9 («i J47r , etc. , 



E 



g(l-g) (Q>+g 2 <c 2 ) 



(i,i;0) 



2/ 2\ ;(2,2;0) 



I N\ 2/ JV\ 2/ JV,JV\ 

g(«i J47T - g («i )4tt + g (« 2 K 



,2/ TV 



2/ W.JVn 



34) 



35) 

36) 
37) 
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The index '47r' denotes quantities that would be measured by the ideal detector with full Air 
acceptance. Substituting into Eq. ( 1761) . we obtain the acceptance scaling formula form Ref. 



VII. FINITE VOLUME CORRECTIONS 

Considering finite system size effects on distributions, we leave the region were the ther- 
modynamic limit approximation is valid. Chemical potentials /ij and temperature T do not 
correspond to the physical ones, which would be found in the GCE, anymore, but have to be 
thought of as Lagrange multipliers, used to maximize the partition function for a given (mi- 
cro)canonical state. First we derive some volume dependent corrections terms, and then find 
a condition that defines the correct values of fij and T. The correct choice allows to write 
down the thermodynamical potentials, the Helmholtz free energy F for CE, and the entropy 
S for the MCE, in terms of the generalized partition function. Some general criterion for the 
validity of the expansion is given. We will compare CE and MCE results with scenarios which 
are accessible to analytical methods in Section IVH CI 

A. Gram-Charlier Expansion 

In Section III II we have shown that in the thermodynamic limit any equilibrium multiplicity 
distribution can be approximated by a Gaussian Eq.(l2TT). Further parameters, describing the 
shape of the distribution, skewness (%), or excess (/C4), tend to zero as volume is increased. 
We return to a generalized version of Eq. (1191) with a number J of conserved quantities and 
use Gram-Charlier expansion |39j. Here we need to explicitly find the inverse square root a -1 
of K2 for the calculation of the normalized cumulants A n Eq. liTTj) . see Appendix [Dl For every 
considered charge one will pick up an additional factor of y/V from the volume element d9j, 



23|: 



acc 



1 — q + quo 



(89) 




00 



—00 




00 



+ 



(90) 



71=3 
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Expanding the exponential in terms of powers in volume, we find: 



Z Q3 



V J / 2 detloi 



n 

3=1. 



2tt 



cxp 



6% 



x 



3! 



V 1 / 2 



+ 



^31,32,33,34 ;\ 

4T 



^"jJ'j/'.rf'j; , 1 - V 



+ 



31,32,33 ^34,35,36 ^-6/ 



2! 3! 



3! 



+ (V- 3 / 2 ) 



(91) 



Correction terms in Eq. fl9Tl) can be obtained by differentiation of exp [— i£-*9j] with respect to 
£ 3 . One can thus reverse the order by first integrating and then again differentiating. Using 
generalized Hermite polynomials, 



= (-!)" ^p 



with the adjusted shorthand notation, 



\ Jl ,J2 J3 



d" 



exp 



(92) 
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the partition function for finite volume can be approximated by: 
Z QJ ~ 



e 2 



27r\/) J/2 det|a| 



1 + + -4r^ + -^r + O (V~ 2 ) 



V 



V 



ys/2 



(96) 



Considering the simplest case of only one conserved charge, it is evident from Eq. (l92|) . that the 
first order correction term in Eq. fl96l) is a polynomial of order 3 in £, while the second order 
correction term is a polynomial of order 4, etc. Hence for large values of £, e.g. a multiplicity 
state far from the peak of the distribution will lead to a bad approximation, and even to negative 
values for P(N). The validity of this approximation is thus restricted to the central region of 
the distribution. We will compare CE and MCE results with scenarios which are accessible 
to analytical methods in Section I VII CI In order to distinguish approximations which include 
corrections up to different orders in volume in Eq. (l96|) . we denote the asymptotic solution 
as CLT (central limit theorem), including terms up to O (V -1 / 2 ) as GC3 (Gram-Charlier 3), 
including terms up to O (V^ 1 ) as GC4, and including terms up to O (V~ 3//2 ) as GC5. 
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B. Chemical and Thermal Equilibrium 



Here we address the question of how to choose the optimal values for T and \x r If exact 
solutions were available, distributions, however not thermodynamic potentials, would be inde- 
pendent of this choice (see as well Section [VHIj) . Our postulate is that our (micro) canonical 
equilibrium state should be as well the most likely state in the GCE. For an isotropic momen- 
tum distribution, the macroscopic state P = is always the most probable state, since all odd 
cumulants involving only momenta vanish n\ x = ^ x ' Px ' Px = . . . = 0. On the other hand, we 
know that the expansion works best around the peak of the distribution. So we choose T and 
fij such that we maximize the partition function at some point equilibrium Q e9 . Taking terms 
up to 0(y~ x l 2 ) into account, the first derivative of the partition function Eq. (IT2"|) reads: 

L („-') \ + 1^ * Hl ( S ) hMM + O (v-) 

"(97) 

The chemical potentials should be chosen such that the first derivative Eq.( l971) of Z® 3 with 



dZ® 3 

dQJ ~ (2 7 r) J/ V( J + 1 )/ 2 det|a| 



respect to the conserved quantities Qi vanishes, hence we maximize Eq.fl96l) at the point Q 3 



eg-. 



. (9£ 



Using only the asymptotic solution, valid in the thermodynamic limit, this condition leads to: 

£ fc = (Qi - V Kl)j ) {a^) j k = . (99) 

Hence the partition function is maximal at the point Q e9 = Vn^j. Charge and energy density 
correspond thus to the GCE values, and /z- 7 — > fi J gce and T — > T gce . While when taking the first 
finite volume correction term in Eq. (1971) into account we obtain: 

-1\ fc 4 



^ (O * + -j^pr ("1, ^ (0 kukM4 = o , (ioo) 

rather than Eq. (l9"9l) . and ^ ^ fi gce , and T ^ T gce . The recipe for calculation of distributions 
for finite volume system thus goes as follows. One should find chemical potentials that satisfy 
condition fl98l) . Then keep in mind that they are chemical potentials only in the thermodynamic 
limit, while for finite volume they are simply Lagrange multipliers. Then one should calculate 
the normalization Z® 3 ' eq and the distribution Z® 3 ' eq using chemical potentials and temperature 
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obtained from Eq. (l98l) . Their ratio gives the distribution P(N\Qi ,eq ) of particles of the selected 
species. A technical comment is in order. From Eq.( l971) it is evident that the first order 
correction term to the derivative of the partition function is a polynomial of order 4 in £, while 
the second one is of order 5, etc. It is therefore crucial to find in numerical calculations the 
correct maximum. 

C. Quality of Approximation 

To test the quality of our approximation for multiplicity distributions at finite volume, 
Eq.( l96i) . for (very) small systems, we compare to analytical solutions for a CE classical particle- 
anti-particle gas, and a classical MCE (without momentum conservation) ultra-relativistic gas. 
The exact solutions are given by Eq. (11211) and Eq. (11361) respectively. Figure [1] shows on the 
top row the multiplicity distribution of positively charges particles for various system sizes in 
the exact form Eq. (I12ip and in different orders of approximation Eq. (I96I) . On the bottom row 
the ratio of approximation to exact solution is taken. In figure [2] the same physical system is 
shown for a (relatively large) positive net-charge. Due to a one-to-one correspondence between 
the distributions of negatively (suppressed) and positively (enhanced) particles we find the 
distribution P(N + ) generally more narrow than in the case of a neutral system. In particular 
towards the edge of the body of the distribution the approximation is worse. For the MCE 
massless gas we compare again approximations to P{N) and ratios to the exact solution on top 
and bottom row of figures [3] respectively for different system sizes. A few general comments 
attempt to summarize these figures. Our first observation is that indeed as system size is 
increased we find a better description of the central region in terms of the asymptotic solution 
CLT (Gaussian with width given by Eq. (l28l) ). Our second observation is that even for systems 
with a very small number (in the order of 5) of produced particles we find a good approximation 
in terms of Gram-Charlier expansion. In particular GC5 provides a very accurate description of 
the central region with deviations in the order of a few percent. This is quite remarkable given 
the fact that multiplicity distributions for such small systems are not smooth and continuous 
functions of multiplicity, while our approximation Eq.( l96l) is. Furthermore we have implicitly 
introduced the concepts of chemical potential and temperature for systems with small particle 
number, which may be in contradiction to the common believe that these parameters can 
only be meaningful when the number of involved particles becomes very large, i.e. in the 
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FIG. 1: top: Boltzmann CE tt+ multiplicity distribution for T = IGOMeV, Q = for V = 100/m 3 
(left), and V = 500/m 3 (right). Exact solution (solid), and in CLT (dash), GC3 (dash-dot), GC4 
(dash-dot-dot), and GC5 (dot) approximations, bottom: same as top, but ratio of exact solution to 
approximation . 

thermodynamic limit. Our last observation is that indeed (see bottom rows of figures [Ul3l) 
finite volume corrections given in terms of polynomials lead only to good results for the central 
region of the distribution. 

To give an estimate for a region in which the approximation is reliable, we note that our 
finite volume approximation scheme begins to break down when the first expansion term in Eq. 
( 1961) becomes unity. In the one-dimensional case this would be: 

^3 (^max 



V 



0(1) 



(101) 
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FIG. 2: same as[TJ but with Q = 5 (left), and Q = 25 (right) 
Approximating the Hermite polynomial if 3 ~ £ 3 , one can get an estimate for £, ma x- 

f o| \ 1/3 

i max ~ - V 1 '* . (102) 



While, when switching back to the definition of £ = > Eq.( fl8|) . the width of the central 

region can be estimated by: 

O \A 3 / 

Hence the width of the central region scales as V 2 ^ 3 , while the width of the curve should scale 
as V x l 2 and our approximation should be quite good. Even though larger volumes work better, 
they will still be sufficiently small enough to allow for calculation of distributions relevant for 
heavy ion collisions. However we want to stress that there is no simple criterion for what 
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FIG. 3: top: Particle number distribution in the MCE without momentum conservation for a ultra- 
relativistic Boltzmann gas. E = lGeV, V = 50/m 3 (left), and E = 3GeV, V = 150/m 3 (right). 
Exact solution (dash), and in CLT (dash-dot), GC3 (dash-dot), GC4 (dash-dot-dot), and GC5 (dot) 
approximations, bottom: same as top, but ratio of exact solution to approximation. 

is a 'small' or a 'large' volume for a particular physical system. Formally the existence and 
finiteness of the (at least) first three cumulants k is sufficient for application of the asymptotic 



solution 39(. Considering the simple case of a multiplicity distribution of BE particles in the 
GCE we find from Eq.()8]) < k^ ,N < k^' N,N < ■ ■ ■ . Hence we find, in particular when finite 
chemical potentials are involved, cumulants growing with order, implying that apart from mean 
and variance further parameters like skewness and excess 38[ of the distribution will remain 
important quantities. 
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D. Temperature and Chemical Potentials in MCE and CE 



The introduction of chemical potentials in the CE and temperature in the MCE was first 
and foremost a mathematical trick which allowed to conveniently integrate partition functions 
for which otherwise no analytical solution could be obtained. However our generalized partition 
function is self-consistent and not in contradiction to the common definition of temperature 
and chemical potential. We will show in the following that our definition of T and fi through 



Eq.( l98l) coincides with expressions well known from textbooks |33f ] 



1. Canonical Ensemble 

The canonical partition function known from textbooks and our generalized version are 
connected as follows (see Appendix IB"]) : 

Z CE (Q,V,T) = ZQ(Q,V,T,^ Q ) e~ Q ^ . (104) 

The Helmholtz free energy F is the thermodynamic potential relevant for CE, 

F = -T\nZ CE . (105) 

Using the first law of thermodynamics dE = TdS — PdV + f^QdQ, and F = E — TS, where P 
is the pressure, and E and S are total energy and entropy, respectively, we can write for the 
differential dF = —SdT — PdV + fiQdQ. The effective chemical potential fiQ associated with 
conserved (electric) charge Q is thus defined by: 

fdF\ # - ^ e-Q^ 

feL S,e-«* ="«■ <106) 



where we have used condition (198]) . Thus the correct choice for the effective chemical potential 
is = 0, which coincides with fiQ = • In the thermodynamic limit this is equivalent 

to, Eq. (l99l) . fiQ — > fiQ.gce- The subscript in Eq. (11061) is used to indicate that the derivative 
with respect to the conserved charge has to be taken at fixed values of V and T. In the non- 
relativistic case, where particle number N, rather than quantum numbers, is conserved, the 
corresponding relation to Eq. (11061) would be (§jj) vt = Hn [33|- F° r a detailed presentation of 
CE calculations please see Section IVIIII 
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2. Microcanonical Ensemble 



The common MCE partition function can be obtain by integration of the GGPF Z ((f) e) and 
multiplication by the inverse Boltzmann factor (see Appendix [C] for details), 

Z MCE (E,V) = Z E (E,V,T) e% . (107) 

The relevant thermodynamic potential in the MCE is the entropy S, 

S ee \nZ MCE . (108) 



The effective temperature is [33fl : 



dS\ + I^el 1 



0^ " T 



(109) 



where we used condition (1981 . = 0. The subscript in Eq. (11091) is used to indicate that the 
derivative with respect to E has to be taken at fixed volume V. Thus Eq. (I109[) resembles the 
optimal choice of an effective temperature for our approximation scheme. In the thermodynamic 
limit, V — > oo, we find T — > T gce , due to Eq.(l99l). 

3. Grand Canonical Ensemble 

Conventionally, e.g. in textbooks, first the MCE is introduced. Summation over energy, 
with temperature being a Lagrange multiplier, used to maximize the entropy, introduces the 
CE. Additionally dropping the constraint of exact charge conservation leads to the GCE. Here 
the chemical potential fig is the Lagrange multiplier. For the MCE (without momentum con- 
servation) and one conserved charge Q the GCE and CE partition functions are defined by: 

Z GCE (V,T,fi Q ) = ^ e Q ^ Z CE (V,T,Q) = £ e Q ^ e~$ Z MCE (V,E,Q) . (110) 

Q Q,E 

While in our notation this line would read: 

z(v,T,fi Q ) = j2 zQ ( v ' t ^q) = E z BiQ (y,T,iiQ) ■ an) 

Q Q,E 

The thermodynamic potential for the GCE is the grand potential Q: 

Q ee -T\nZ GCE = -T\nZ\^ =0 . (112) 
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Again, if exact solutions of the canonical or microcanonical partition functions were available 
this reversal would not have been necessary. However this redefinition of the GCE partition 
function is entirely consistent and simplifies calculations considerably. Whenever an exact 
solution to our generalized partition function is possible, all the above relations would hold 
exactly. We believe this interpretation of the GCE partition function as the generating (or 
characteristic) function of a statistical system to be quite useful, even in more general cases 
than the one presented here. 



VIII. THE SIMPLEST EXAMPLE 



A. An Exact Solution 



To clarify the procedure described in the previous sections, let us consider a simple example 
of an ideal Boltzmann particle-anti-particle gas (q+ = 1, g_ = —1). In this simple case the 
calculations discussed above can be done explicitly. The GCE partition function reads: 



n_|_ ri- 



Z = —\—\ = Yl z ( n +> n -) = ex P [ z + + z - ] = ex P 2z cosh (^fr) 

71j . L \ _Z / . 

n+ ,n- =0 n + ,n- =0 

(113) 

where z± — z exp(±/j,q/T), and z is a single particle partition function in MB approximation: 

gV 7 2 / y/p 2 + m 2 \ gV 2 (m\ 
z=—J P dpe W \ I = —mTK 2 { ¥ ), (114) 

o ^ ' 

where g and m are respectively the degeneracy factor and particle mass, and K 2 is the modified 
Hankel function. The conserved charge Q is just the difference of N + and AL. To be definite 
we discuss the distribution P(N + ) of positively charged particles. In the GCE it can be easily 
found: 

all states with iV + particles Er + ,n_=o z ( n +> n -) 6 ( N + ~ n +) 
{ +j all states ~ E~ + ,„_=o^( n +.«-) 

N + 

= exp(-^ + ) . (115) 

As can be seen it has the form of the Poisson distribution with the following first two moments: 

oo oo 

(N + ) gM , = ^iV + P(iV + ) = z+, (N 2 + ) g . c . e . = J2 N l P ( N +) =*+ + *+■ (H6) 

N + N+ 
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This gives the GCE scaled variance: 



CO 



(Nl) 9 . c . e . ~ (N + ) 



g.c.e. 



g.c.e. 



1 



• ' (117) 

In the thermodynamic limit z — > oo the Poisson distribution (11151) can be transformed into the 
Gauss one. Using Stirling's formula, N + \ ~ y/2irN + exp (N + \nN + — N+), for N + 1 one 
finds at \N + — z\ -C z: 



P(iV 4 



. AT 4 



NJ 



1 



exp(— z^ 



exp 



2z, 



118) 



V 27r z -\ 

The distribution of net charge P(Q) and the joint distribution P(N + , Q) can be found by the 
straightforward calculations: 



1 A _ r „ zl + z n _~ e«^/ T 



P(Q) = y E W-K-n-]) 



n+ ,n_=0 

oo 



n i ! n_! 



P{N+,Q) = - £ *(Q-[n+-n-]) 



n+,n_=0 



71+ 71- 



n+! n_ 



(119) 



2JV+-Q 



Z iV + !(iV + -g)! ' 

(120) 

where Iq(2z) is the modified Bessel function. The Eqs. (11191 1120| are the simplest appearance 
of the general Eqs. (TSJ HJ). One can also notice that for this example the CE partition function 
from the Eq. (T5J) equals the modified Bessel function = Iq(2z). 

From the Eqs. (11191 11201) one finds the CE particle number distribution (see the Eq. ([HJ)): 



P(N+\Q) 



P(N+,Q) z N + z N +- 



[i Q wr 1 



P(Q) N + \ (N+-Q)\ 
The first and second moments of the CE multiplicity distribution can be easily found [28]: 

£ N + p ^ = * w ■ 

N+=Q yv ; 

f; ni p(n + \q) = z + <- 2 Iq - 2{2z 

N+=Q 

This leads to the CE scaled variance: 



121) 



(N + ) 



2 

c.e. 



<#+>c 



1 - z 



I Q (2z) ' I Q (2z) 
Iq-x&z) I Q - 2 (2z) 



(122) 
(123) 



In the thermodynamic limit when z — > oo and Q/2z = y one finds 

1 

: 2 ~ ~2y/T+tf 



z(y + v 7 ! + y 2 ) 



y 



(124) 



(125) 
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Let us compare the results for the average particle number (N + ) and for the scaled variance 
lo + calculated in the GCE and CE. We fix the V and T parameters to be the same in the GCE 
and CE. The relation between the chemical potential \xq of the GCE and the charge Q of the 
CE formulation is obtained from the following requirement: 

Q = (Q) = (N + ) g . c . e . - (iV_),. c . e . = z + - z_ = 1z sinh (^) . (126) 

This gives, exp (/xq/T) = y + y/i + y 2 , and leads to (iV + ) c . e , ~ (vV + ) a . c . e ., which means the 
thermodynamic equivalence of the CE and GCE. Comparing Eqs. ( I117P and (j!25p one finds 
00 t.e. u g.c.e.i thus, the scaled variances are not equivalent even in the thermodynamic limit. 

We have used the GCE distributions P(Q) and P(N + ,Q) to calculate the CE distribution 
P(N + \Q). The Eq. (I12ip demonstrates that the /iQ-dependence has completely disappeared in 
P(N+\Q) if exact analytical results in the GCE for P(Q) (11191) and P(N + ,Q) (11201) are used. 
This means that the GCE with arbitrary value of /xq (and, thus, the corresponding value of 
(Q)) can be used for the exact CE calculations of P(N + \Q). 



B. The Saddle Point Expansion 



We will discuss now the approximate calculations of P(Q) (11191) and P(N + , Q) (I120p in 
the thermodynamic limit. We replace the corresponding delta- functions which fix the value of 
net charge, Q, and positively charged particle number, N + by their Fourier representations as 
the 4>q and <p + integrations. This method can be used for more complicated cases when exact 
analytical results can not be obtained. In these asymptotic calculations the role of the chemical 
potential will be demonstrated. The Eq. (11191) can be rewritten as follows: 



Z 



(127) 



1 

~Z 



2ty 



exp( 



exp 



z exp + z<; 



+ z exp 



T 



In thermodynamic limit, z — > oo, one can expand exp (± [/xq/T + i(j)Q[) in the Taylor series 
and leave only the terms up to <f>Q, because for z — > oo the main contribution comes from the 
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= region. Then the distribution P(Q) becomes a Gaussian: 

oo 

P(Q) = j y^ exp [-i (Q - 2z sinh(/i Q /T)) Q - z cosh(/i Q /T) + 



~ P G (Q) = [47rzcosh(/i Q /T)]- 1/2 exp 
For P(N + , Q) at z — > oo one similarly finds: 



(Q - 2z sinh(/x Q /T)) s 
4^ cosh(/i Q /T) 



(128) 



P ^ = Z ./ 2. 



2tt 



TT 7T 



Z 7 2n J 2tt 

— 7T — 7T 



exp 



—i 



- iN + <p + + 2 exp ^ + 20 Q + + zexp (-^ - ^ 



oo oo 

Pc(iV+ ' g) = z J 17 J 17 exp 



- 0+ + y - z + y - z + 



(129) 



The integration over 0q and 0+ in Eq. (11291) gives: 



Pg(N + ,Q) 



2ixz 



exp 



iexpfo/T) (Q - {Q)f + -exp(// Q /T) (Q - <Q» {N+ - (N+)) 



-~cosh(// Q /T) (iV + -(iV + )) 2 

z 



where (Q) and (iV + ) in Eq. (I130j) correspond to the GCE values, 

(Q) = 2z sinh(/i Q /T) , (N+) = z exp^/T) . 



(130) 



:i3i: 



The Eq. (11301) has the form of Bivariate Normal Distribution i.e. Gauss (Normal) distribution 
in two dimensions [38]. The CE distribution (11211) is then approximated as, 

P(N + \Q) = - Po(N + \Q) S 



p{Q) 



+ 



Ti z cosh 1 \pq/T] 
exp(/x Q /T) 



exp 



Pg{Q) 

cosh (un/T) ., r ,„ r . .o 
; (JV+ - (N + )) 2 



(N + -{N + ))(Q-(Q)) 



exp(WT) (Q-(Q» 2 



2^: cosh (//q/T) 



(132) 
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Comparing the distributions (11281) and (1130)1 with the exact expressions for P(Q) ( 11191) . and 
P(N + ,Q) ( 11201) one can notice that the Eqs. f l 11 9 j) . and ( 1120p contain the dependance on the 
chemical potential hq just as the factor e^ Q ^ T . Thus if one succeeds in the exact calculations 
the resulting CE distribution P(N + \Q) f 1 1 2 1 j) does not include /xq dependence in contrast to 
(11321) . It means that the choice of the chemical potential is irrelevant for the exact CE calcu- 
lations (usually it is chosen equal to zero), while the value of hq is crucial for the approximate 
calculations. The saddle point expansion works if \iq chosen to fix (Q) = Q in thermodynamic 
limit, i.e. for Hq/T = arcsinh(Q/2z). Thus the GCE should be thermodynamically equivalent 
to the CE with fixed Q net charge. 

Let us illustrate these statements. We plotted the exact (11191) and approximate (Gauss) 
distributions P(Q) (11281) for the arbitrarily chosen values z = 20 and Q = 50 with zero chemical 
potential hq = (Fig. HI left) and for the chemical potential hq ^ that corresponds to the 
condition (Q) = Q (Fig. left). 
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FIG. 4: The exact and approximate distribution P{Q) (left), their ratio (center), and the integrant 
(Int) from the Eq. (|128j) (right) for the 'wrong' chemical potential ^q/T ^ arcsinh(<5/2z). 
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FIG. 5: The same as in Fig. H]but for correct fiQ = arc sinh(Q/2z). 



One can see that if one chooses the wrong chemical potential /ig = than the approximate 
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distribution (11281) is approximately 30 times smaller than the exact one (11191) at the point 
Q = 50 (see the ratios Pg(Q) I ' P{Q) Fig- HI El center). For non zero chemical potential fiQ ^ 
that satisfy the condition (Q) = Q, i.e. for the chemical potential Hq/T = arcsinh(Q/2z), 
the exact and the approximate distribution almost coincide in a wide range near Q = 50 (see 
Fig. [5], left and center). 

The mathematical meaning of this physical requirement, (Q) = Q, is the follow- 
ing. The approximate formulas always have an oscillating part exp [iA(p — B(f> 2 ] d(p = 
cos [A(f\ exp [— B(p 2 ] d(p (the imaginary part includes sin[A0] that gives zero after in- 
tegration because of symmetric integration bounds). This oscillating part is nonzero 
4 / for (Q) 7^ Q thus we need to take the terms higher than 2 in or- 
der to make a good approximation. The integrant (Int) from the Eq. ( 11281) . Int = 
exp [— i (Q — 2z sinh(/iQ/T)) 0q — z cosh(/iQ/T) 0q + ...], is shown in the Figs. HI [5] (right). 
One can see that for (Q) ^ Q (Figs. HI El center) the Int oscillates around zero fast, thus even 
small mistake in numerical integration, that is unavoidable in real calculations, gives the wrong 
result for P(Q). However if we take the correct fig that satisfy the condition (Q) = Q, then 
Int = exp [— z cosh(yUg/T) 0q + . . .] becomes smooth and the term 4>q is enough (see Figj5l 
right). 

At the point (Q) = Q from the Eq. (11301) one finds: 



Pg(N + \Q) 



ttz cosh 1 (/Xq/T) 

1 



exp 



(N + ) 



exp 



(JV + - (iV + )) 2 
z cosh -1 (hq/T) 

(N + - (N + )) 2 
2 w+. • (N+) 



(133) 



In Eq. (11331) the asymptotic values for multiplicity (valid in the GCE and CE), (N + ) 



zexp [jmq/T], and for the scaled variance in the CE, = [2 exp (/ig/T) cosh (/ig/T)] -1 (see 
Eq. (11251) ). have been used. 

The above procedure can be generalized for the case of n integrations using the following 
formula: 



n 



2tt 



exp 



iA T 4> 



1 



B 



exp 



-\ A T B- l A 



(27r) n /Vdet| J B| 



(134) 



where T means transposed vector and B 1 is the inverse of a nonsingular matrix B. For 
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example, the distribution P(N + , Q) has the following parameters: 

A = h~ Z -- Q \, <?=M, B=h + Z - Z+ ), Aet\B\=z + Z-=z\ (135) 

V -+~ N + J \<f>+) V 2 + z +/ 

One can easily check that the substitution of the Eqs. (11351) into Eq. ( 11341) gives P(N + , Q) from 
the Eq. (TT30D . 



C. Microcanonical Ensemble 

re same way. In the simplest 



The microcanonical partition function can be calculated in t 
case of massless Boltzmann particles we know the exact answer 



27|]: 



1 1 ' " P(E) " Z MCE E(3N -1)\N\ ' 1 J 

where Z MCE = ^ 0-P3 (; |, §> 2; J^) is the MCE partition function, F 3 is a generalized hyper- 
geometric function, and x = gVE 3 /ii 2 . Let us consider P(E) in Eq. (I136p . One can obtain it 



by straightforward iterative calculations similar to the Ref. 



271]: 



iv=i • v ; ./ \ ) j fc=1 fc=1 

N=l • \ / JO J k=l 



e/t 7 MCE 

' Z - (137) 



Z 

The saddle point expansion for P(E) gives 



7-oo 2vr £oo2ttT 



Fc(£) = 7™ e ^ ' (138) 

where (N) = gVT 3 /Tr 2 , (E) = 3(N)T and Z = e {N) . The exact distribution P(E) and its Gauss 
approximation P G {E), their ratio R = P G (E)/P(E), and the integrant Int = R e [e^ N ^- E / T ^ E \ 
(the imaginary part vanish because of symmetric integral bounds) are shown in the Figs. [6], [7| for 
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the parameters E = 20 GeV, V ~ 768/m 3 and for two different temperatures: the 'equilibrium' 
temperature T = 160 MeV that gives (E) = E and for 'non-equilibrium' temperature T = 
120 MeV that gives (E) = 6.32 GeV^ E. 




E (GeV) E (GeV) 4> E < rad > 



FIG. 6: The comparison of the exact distribution P(E) and its Gauss approximation Pq{E) (left), 
their ratio R = Pq{E)/P{E) (center) and the integrant Int = Re[e^ 3( - N ' ) ~ E ^ T ^ 1 '' fiE ] for the parameters 
(E) = 20GeV, T = 120 MeV 




E (GeV) E (GeV) <I> E (™« 

FIG. 7: The same as Fig. [6] but for T = 160 MeV. 

Comparing the Figs. HI [5] and Figs. El [7] one can see that if one choose the wrong value 
of temperature than the integrant Int in P(E) oscillate fast and Gauss approximation gives 
wrong result for the MCE system see Fig. [6J It means that if we calculate the MCE distribution 
P(N\E) by means of the auxiliary GCE distributions P(N, E) and P(E) than the temperature 
play the role of the 'chemical potential' similarly to the CE. 

D. Resonance Decay 



In this subsection we would like to provide a connection between probability distribution P, 
cumulant generation function and the generating function G 26j that was firstly used to 
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take into account resonance decays analytically. Let us consider again a simple example of an 
ideal Boltzmann gas of positively and negatively charged particles, i.e. particle-antiparticle gas, 
and, additionally, 'allow' the decays of positively charged particles as well as their antiparticles. 
The formula for P(Q), (I127p . remains the same, while P(N + ,Q), (11291) . changes: 



- iN +( p + 



7T 7T 

u/at n\ 1 [ d ^Q [ d< ^+ 

p(iV +' g) = zj ^rj ^r exp 

+ z exp (Jy + i( Po) Yl r + exp c< ^+) + z exp ir^T ~ i( ^ Q ) 5^ 6XP ^ ° q 



(139) 



where c is a multiplicity of a particular decay channel and is a branching ratio that corre- 
sponds to the decay of positively/negatively charged resonance. 

There is no any real system that corresponds to the partition function (11391) . but it is the 
simplest example where we can show how these formulas work. 

We can expand the exponent in the Eq. (1 139ft : 

P(N + ,Q) = ^j^Rj^± exp[-*(Q-[n + -n_])0 o - ziV + + ] 



x 



E 



2. Z_ 



exp u cq 



n+,n_=0 

and immediately recognize the generating function that was used in our previous work 



(140) 



26): 



T+ exp (i ccf> + ) r - ex P (» C( f 



The difference from the generation function used in [26( is that here we study the case with 
only two types of resonances, n±, and we also use the different definition for a decay channel. 

From the Eq. (11401) one can also see that the cumulant generating function ^ is connected 
with the generating function G in a very simple way: 



In 



E 

n_i_,n_=0 



n+ n- 
Z, Z_ 



nil ri-l 



G(0 H 



(142) 



In order to make a comparison with the Ref. [26j one should also remember some technical 
moment that in the paper 26( we used the auxiliary parameters A that are the same as exp(i(j) + ) 
in this paper, however 







A- 



_d_ 



and 







A, 



_d_ 



x. 



_d_ 



(143) 
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From the Eq. ( 1142ft it also can be easily seen that we can multiply the generation function by 
any other normalized distribution for example by the finite acceptance probability distribution 
and it will enter the partition function in the same way as the generating function G (see 
Section EE}. 

IX. SUMMARY AND CONCLUSION 

We have presented an analytical expansion method for calculation of distributions at finite 
volume for the canonical as well as the microcanonical ensembles of the ideal relativistic hadron 
resonance gas. The introduction of temperature into the microcanonical partition function and 
chemical potentials into the canonical partition function have lead to the identification of the 
grand canonical partition function with the characteristic function of associated joint probability 
distributions. The microcanonical and canonical multiplicity distributions could thus be defined 
through the conditional probability distributions of finding a certain multiplicity while other 
parameters (global charge or energy) were taken to be fixed. 

We have shown that in the thermodynamic limit the central region of the canonical and 
microcanonical partition functions can be approximated by multivariate normal distributions. 
Multiplicity distributions tend thus to Gaussians as system size is increased. In particular, we 
could find a simple formula for the scaled variance, hence the width of this Gaussian, in this 
limit. It was further possible to show the equivalence of our results to those obtained from the 
microscopic correlator approach. 

Apart from providing a simplified formula for asymptotic multiplicity fluctuations, the ap- 
proach presented in this work owns a few conceptual advantages. In considering finite volume 
corrections to the system partition function, thus relaxing the assumption of thermodynamic 
equivalence of different statistical ensembles, one is lead to demanding that the partition func- 
tion should be maximized for a particular set of conserved charges. It turned out that this 
requirement is entirely equivalent to the well known textbook definitions of chemical potential 
in the canonical ensemble as the derivative of Helmholtz free energy with respect to conserved 
charge and temperature in the microcanonical ensemble through differentiation of entropy with 
respect to conserved energy. 

Our method is based on Fourier analysis of the grand canonical partition function. Conven- 
tionally one would not introduce chemical potentials and temperature into these calculations. 
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However one then faces the problem of a heavily oscillating (or even irregular) integrant, making 
numerical integration unpractical. Artificially introduced temperature and chemical potentials, 
correctly chosen, produce a very smooth integrant allowing for expansion of the integrant in 
powers of volume. Analytical solutions to asymptotic multiplicity distributions could thus be 
found in terms of Laplace's expansion, while finite volume corrections could be obtained from 
Gram-Charlier expansion. A first comparison with available analytical solutions to simple sta- 
tistical systems suggests that good results can be expected even for rather small volume. One 
drawback is that the results can only be applied to the central region of the distribution, owning 
to the fact that finite volume correction terms appear in the form of Hermite polynomial of low 
order. 

Additionally we have included resonance decay directly into the system partition function 
of the hadron resonance gas model. This treatment has proven to be more economical than 
the previous handling through generating functions or Monte Carlo techniques. Neglecting 
correlation in momentum space a first order approximation to the effect of finite acceptance 
can be made. 

This paper provides a connection between the two previously published methods for cal- 
culation of multiplicity distributions, the microscopic correlator approach and saddle point 
expansion, and should be of great practical use. 
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APPENDIX A: CHARACTERISTIC FUNCTION 



Let us consider a probability distribution of some observable \ which can take values f x (x) 
on the real x axis with the following normalization: 



dx fx( x ) = L 



(Al) 



Hence f x (x) defines the probability of finding a value x for the observable x- One can define the 
characteristic function <$> x (iO) as the Fourier back transformation of the probability distribution 

fx( x ) : 



dx e lx0 f x (x) . 



(A2) 



The moments £ n {x) can directly be calculated from the characteristic function from the re- 
spective derivatives at the origin: 



d n $ x (U 



dfr 



8=0 



i n / dx x n f x {x) . 



(A3) 



Once the characteristic function & x (i9) is known, the probability distribution f x (x) can be 
found by Fourier transformation: 



7T 



(A4) 



APPENDIX B: CE PARTITION FUNCTION 



The (generalized) system partition function of a system with conserved charge vector Qi 
(B,S,Q), Eq.©, reads: 



z QJ 



e~ iQJ ^ exp 



E 



z i WjifJ-j, 



The single particle partition function Eq. (jSJ) includes chemical potentials as well as Wick rotated 
fugacities. Expanding the logarithm in Eq.® yields: 



Zl 



9iV 
(2vr) 3 



dp In (lie t e t e 



±1 



(B2) 



(B3) 
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where the effective chemical potential of particle species I is given by \i\ = erf Hi, and the 
single particle partition function of a 'lump' of mass n\ mi of species / is given by 32j z ni = 
-B^rg J d 3 p n^ 1 e^'r 1 . Eq. flBlj) can be thus written as: 



z Q3 



3 

— 7!- 



e -iQH 3 exp 



l n;=l 



nnz 



((=F1)"' +1 2 „,)' 



\ni+l \»nj 



(B4) 
(B5) 



re; • 



Z n;=lfc„ ; =0 

Only sets of numbers that meet the requirement: 

oo 
i nj=l 

have a non-vanishing contribution to the integral. For any such set one finds: 



(B6) 



i n,=l 



(B7) 



Therefore one can pull the factor Eq. (1B7I) in front of the integral flBlj) : 



Hit 



-iQi<t>. 



00 00 (~ \k ni 

Yin j2 c ikn ' niq3A < 



(B8) 



However this is just the standard definition of the system partition function times the factor 
Eq.jEZJ, thus: 

J oo 



z Q3 



and 



jGCE 



n e 

j=l QJ=-oo 



(B9) 



A similar calculation can be shown to hold for the MCE. 



APPENDIX C: MCE PARTITION FUNCTION 



Having introduced temperature in the MCE partition is probably a somewhat unusual ap- 
proach. Here we consider this in detail in the example of a massless gas in Boltzmann approxi- 
mation without momentum conservation. We will first state and solve our version of the MCE 



partition function and compare to the analytical solution known from textbooks (e.g. 331 ]) 
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We adopt the notation of [27], Z E ' N = W N (E, V) and Z E = W (E, V). The number of states 
consistent with the constraints of fixed energy E and particle number N in a GCE system is 
given by the double Fourier integral over its GGPF Z {4>Ni 4>e)'- 



W N (E,V,T) 



2tt 



2tt 



- e~ iN(t>N e~ tE4>E exp 



— / p 2 dp e" ft ' ' "'" L ' ' v ' x 



(CI) 



Expanding the exponential and solving the integral over 0tv yields: 



W N (E,V,T) = -^ f ' 



N\ V2vr 2 



N 



2tt 



-i N 



p 2 dp e 



(C2) 



Further now solving the integral over momentum gives: 



W N (E,V,T) 



i fgv 



N\ V27r 2 



N 



>E_ e -iE<t> E 



2tt 



i N 



(C3) 



Eq. ( 1C3I) has obviously a pole of order 3iV at <pE = —i/3. So we close the integration over the 
lower hemisphere and use the residue theorem. 



-i fgV 



W N (E,V,T) = — Res 



it 



N 



-iEcj> E 



3N ' 



<p E = -i(3 



-%_ fgV 



it 



N 



lim 



gv 



71 z 



N 



(3iV-l)! 4>E^-m 
E3N -' .-/w 



d 3N - 



- (<Pe + 



3N 



-iEcj>L 



{(3 - ifa 



,3V 



AT! (3JV- 1)! 



(C4) 



M ultip lication with the inverse Boltzmann factor yields the result known from the literature 



Vlultipnc; 



331, W N (E, V) = W N (E, V, T) e? E : 



W N (E,V) 



gv 



N 



E 



3JV-1 



AT! (3iV- 1)! 



(C5) 



Using condition Eq. (|98[) it further follows that the temperature is given by the average energy 
per particle T = E/(3N — 1). In complete analogy to what was presented before: 

oo oo 

W(E,V) = Y^W N (E,V), and W(E,V,T) = ^ WW [E, V, T) . (C6) 

7V=1 N=l 

However in both case we define the MCE multiplicity distribution by: 



P(N\E) = W N (E, V) W~ l (E, V) = W N (E,V,T) y\T x (E,V,T) 



(C7) 
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The introduction of temperature in Eq.f lCll) does not necessarily simplify this calculation (it 
even drops out entirely from the r.h.s of Eq.( 1C7l) ). however makes approximations possible. 
The GCE partition function finally is given by: 



oo oo 

Z GCE {V,T) = 1 + J dE W(E,V)e- /3E = 1 + J dE W (E, V, T) = exp 



° VTS ~ .(« 



7T 2 





The additional term ' + 1' arises form the fact that the MCE partition function only includes 
states with at least one particle, while the state of zero energy and zero particle number con- 
tributes to the GCE partition function. 

APPENDIX D: WIDTH AND NORMALIZATION 

The inverse of a nonsingular J x J matrix R<i can be obtained from its adjoint and its 
corresponding minors |38| : 

( ~ n = adjfc, = [ ) T (D1) 

v 1 det|« 2 | det|/c 2 | 

In our case only the minor = det \ k2\ is of interest, hence we find for the element in the 

upper left corner of the inverse of k^. 

(*•-'),. ~ W ^ 

Generally we can find the inverse of a matrix by multiplication of a diagonal matrix, with the 
inverses eigenvalues on its diagonal r~l, with ortho normal transformation matrices, which can 
be formed from the eigenvectors v a of k%, 

= E E < D3 ' 

a=l b=l 

where Vi >a is the a th component of the i th eigenvector, vj^ is the transpose of Vi >a , hence a matrix 
with eigenvectors in its columns, and {T~ 1 ) a b is a matrix with the corresponding inverse of the 
a th eigenvalues t~ l on its diagonal. We find for the upper left most entry: 

(^)m - E < c - ^g- (~) 

a=l 1 Z| 
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Likewise, we can express the inverse sigma tensor in term of eigenvalues and eigenvectors of K2 
(see Section IVlTl) : 

(*-)„ - EE «u i-~ 1/2 U <> = E «u fl») 

a=l 6=1 a=l 

We just need the sum of squares of the entries of the left most column: 

E = E E E ^ ^ *: 1/a ^ «^ *r 1/3 = E < ( D6 ) 

i=l i=l a=l 6=1 o=l 

since ^1=1 v^v^f, = 5 a ^. This coincides with Eq. (1D4j) and proofs Eq. (l27l) . 



i=l 



in the most general case of a J dimensional 2 nd rank tensor K2- 
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